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ABSTRACT 


The Department of the Navy (DON) has prioritized the establishment of a robust 
renewable energy platform for its installations. Solar energy has emerged in a way that 
has provided viable alternatives. The use of solar energy, however, poses a problem of 
maintaining continuity of power during hours when sunlight is not available. The 
development of Thermal Energy Storage (TES), while addressing this problem, 
introduces complexities that can reduce overall plant efficiency. This research aims to 
simplify the plant design in a way that minimizes such complexities and maximizes 
efficiency. 

The research in this thesis explores a new design that integrates static molten 
metal baths with the solar collector. The proposed initial design is centered on a simple 
open air Brayton cycle. Aluminum was selected as the phase change metal, providing 
latent heat to the air in the collector during hours of non-daylight through freezing. The 
research focuses on the thermodynamic analysis of the system, initial collector design 
requirements, and minimum volumetric requirements of the aluminum phase change 
metal. An analytical and numerical flow optimization study of the system was performed 
for comparison to design calculations. 
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I. 


INTRODUCTION 


The purpose of this thesis is to present an analysis of a new solar thermal plant 
design that integrates Thermal Energy Storage (TES) into the solar collector heat 
exchanger. The availability of sunlight is intermittent; therefore, TES systems are 
required to supplement the solar energy in order to meet the demands of supported 
installations. TES as separate sub-systems, however, can introduce complexities that are 
costly and reduce the overall efficiency of the plant. The design of the plant in this thesis 
focuses on simplicity in a way that eliminates the complexities presented by current TES 
systems. The goal is to provide a launch pad for continued research and refinement of the 
design. 

The standard Brayton cycle describes the basic design of the plant. The working 
fluid, air, receives heat directly via the concentration of the sun’s radiation into the solar 
collector heat exchanger. The heat exchanger air ducts in the prototype design are 
embedded in a static molten aluminum bath. This bath requires no pumping of the metal 
and simply freezes and thaws around the air ducts. The latent heat energy from this phase 
change metal continues to heat the air as the metal freezes. This thesis provides a 
thermodynamic analysis of the plant, minimum volumetric requirements of the aluminum 
phase change metal, as well as basic geometric design requirements of the solar collector. 

In addition to the prototype design requirements, the air flow was optimized to 
observe plant behavior. Given a specified design power, a basic thermodynamic analysis 
yields a design mass flow specification. By optimizing the mass flow for power, 
however, the same plant design could potentially exceed design power, allowing for a 
transitional operating range between day and night. Analytical and numerical methods 
were conducted in this study to determine the effects on power while varying mass flow. 

A. ENERGY GOALS 

The United States Department of Defense (DOD) is the largest consumer of 
energy of all the federal government agencies, making up at least 80% of the 
government’s energy consumption, as stated by Schwartz [1]. Of the DOD’s energy use, 
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Schwartz [1] also included that installations make up 25% of that share. In response, the 
DOD has taken action to revise its strategy toward energy, investing heavily in the 
research and development of renewable energy sources used to supplement commercial 
grid power for defense installations. Ultimately, the goal of the DOD, as discussed by 
Feldman [2], is to produce 25% of its energy from renewables by the year 2025. 

The Department of the Navy (DON) has responded to the DOD’s energy goals 
with a more aggressive agenda. In the DON’s Strategy for Renewable Energy, Secretary 
of the Navy (SECNAV) Ray Mabus set five energy goals for the Navy, one that directly 
applies to energy use by installations. He directed that “by 2020, DON will produce at 
least 50% of shore based energy requirements from alternative sources; 50% of DON 
installations will be net-zero” [3]. Security and independence were asserted by Mabus [3] 
as the Navy’s two energy priorities. Based on these priorities, the DON has taken action 
to reduce its carbon foot print both operationally and with its installations. 

The DON currently leads the DOD in solar energy development, as shown by the 
Solar Energy Industries Association (SEIA) [4]. Eigure I illustrates the comparison of the 
Navy to other branches of service with regard to installed solar energy. 



Navy Air Force Army 


Eigure 1. Installed Solar Capacity by Branch. Source: [4]. 
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By comparison, the Navy has made significant gains over its counterparts in the 
development and implementation of solar energy. As of 2013, due largely in part to the 
implementation of the Navy’s One Gigawatt (IGW) program, the Navy has installed over 
58MW for various installations, as stated by SEIA [4]. It is through this program that the 
Navy executes its goals set by SECNAV. According to McGinn [5], by 2015 the Navy 
had actually exceeded its energy goals for installations, but not entirely due to organic 
renewable energy production on its own installations. Part of the IGW program included 
provisions known as Power Purchase Agreements (PPA) that allow the procurement of 
renewable energy from outside commercial sources, according to the DON [6]. While 
PPAs might be a necessity for some installations where land availability is a challenge, 
they do not allow installations to operate independently from commercial grids. Eor true 
energy security, an installation would need to produce its own power through organic 
growth of its power production. 

Installations continue to be a viable candidate for the development of solar 
energy. Unlike vessels and other forms of operational equipment where space is limited, 
shore based installations provide the necessary space to build full-scale plants capable of 
supplying a significant share of the demand. 

While solar power is promising in enabling the DON achieve its strategic goals, it 
does not come without challenges. The primary concern of utilizing solar energy to 
power facilities is overcoming the intermittent loss of the power source, the sun. It is 
during these times that the disruption of power must be addressed, especially as a matter 
of security for DOD installations. 

B. BACKGROUND AND LITERATURE REVIEW 

The sun’s radiation is an excellent source of clean renewable energy. Many 
factors, however, can cause disruptions in the power source such as non-daylight hours, 
adverse weather conditions (cloud cover), and seasonal variations throughout the year. 
These factors are more prominent in some geographical regions than others that can 
further complicate the planning and development of solar renewable energy in DOD 
installations across the nation. 
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The National Renewable Energy Laboratory (NREL) [7] published a map of the 
United States that illustrates by color code the variance of concentrating solar resource on 
an average annual basis. Eigure 2 illustrates how the change in local geographic climate 
across the country can provide either more or less advantage for the utilization of 
concentrating solar thermal energy. 
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Eigure 2. Concentrating Solar Resource of the United States. Source: [7]. 


The southern regions of California and Nevada, as well as Arizona and New 
Mexico, provide the best prospects for the development of solar power due to the 
generally dry and sunny climates. In other regions weather and seasonal variations are 
more prominent that result in a greater need for efficient TES systems to supplement the 
solar energy source. 

Thermal Energy Storage refers to methods where thermal energy is extracted and 

stored for later use, as described by the International Renewable Energy (IRENA) [8]. 

Many design concepts to supplement larger scale concentrating solar power plants 

currently exist that achieve this purpose and typically fall into two main categories: 

sensible heat storage and latent heat storage, both of which were discussed by IRENA 
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[8]. Figure 3 is an example of a centralized solar thermal plant that uses thermal salt 
thermal energy storage. 



Figure 3. Centralized Solar Thermal Plant with Thermal Energy Storage. 

Source: [9]. 


IRENA published general efficiency ratings for TES modes in Table 1. TES 
systems are not 100% efficient. All forms of TES will lose some heat to the surroundings 
that can be mitigated through proper insulation technologies. Most forms of TES, 
however, require energy in the form of pumps or other devices to circulate the energy 
back into the cycle for continued power. Heat losses and design complexities can 
therefore impact the efficiency of these systems. 
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Table 1. Efficiencies of Thermal Energy Storage Systems. Adapted from [8]. 


TES Mode 

Efficiency [%] 

Sensible Heat 

50-90 

Eatent Heat 

75-90 


1. Sensible Heat Storage 

Sensible heat TES refers to methods where energy is transferred to a medium 
without the involvement of a phase change. Many designs using this method have two 
reservoirs of hot and cold media. Common sensible heat storage systems use water or 
molten salts as the heat transfer medium. 

While the cost of such systems may be relatively low, the volumetric 
requirements are higher due to the lower energy density of these materials, as discussed 
by IRENA [8]. These systems require additional energy to circulate the media in order to 
continue powering the cycle that can lead to lower plant efficiency and increased plant 
complexity. 

Herrmann and Kearney [10] show how sensible heat TES systems are among the 
most widely used form, with molten salts being the preferred candidate for larger scale 
plants. Other forms have been investigated, such as water. While hot water storage is 
inexpensive and easy to procure, its use is typically limited to micro grids, as discussed 
by Hinke [11]. Water has a low energy density; therefore, the resulting large volumetric 
requirements could be problematic when scaling up for commercial use. Eor this reason, 
existing commercial solar plants have widely been using molten salts, such as the 
Andasol Power Station in Eigure 4, located in Andalusia, Spain. The Andasol Power 
Station in Spain operates a molten salt sensible heat thermal energy storage system. It 
consists of three power plants, each rated at 50 MW, and provides heating from thermal 
storage for up to eight hours [12]. 
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Figure 4. Andasol Parabolic Trough Solar Power Plant with 
Molten SaltTES. Source: [13]. 


Current designs of molten salt TES systems are complex, involving the circulation 
of up to three heat transfer mediums: a heat transfer fluid (HTE), typically oil, the molten 
salt, and a Rankine cycle working fluid (water), as stated by Herrmann and Kearney [10, 
14]. Eigure 5 provides an illustration of a typical molten salt TES system powering a 
Rankine cycle with a parabolic trough collector. 



Eigure 5. Parabolic Trough Solar Power Plant Utilizing 
Molten SaltTES. Source: [15]. 
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In this system, the HTF is heated through the parabolic troughs that direct the 
sun’s radiation to heat transfer pipes. The hot HTF is circulated to power the Rankine 
cycle as well as heat the molten salts that are circulated through a common heat 
exchanger, shown by Abutayeh et al. [14]. During times when the sunlight is not 
available, the molten salts are then recirculated while heating the HTF to power the 
Rankine cycle, as discussed by Abutayeh et al. [14]. 

Complexity is often the main disadvantage of using molten salts as thermal 
energy storage, as well as the added volumetric requirements, although efforts have been 
made to explore design variations that address these problems. For example, Gil et al. 
[16] proposed a new concept that reduced the molten salt storage volume to one tank that 
is directly heated by the sun. The design prototype would also remove the need for 
pumping the salts as a result. 

2. Latent Heat Storage 

Other forms of TES systems utilize the latent heat of phase change materials 
(PCM) to provide additional thermal energy needed to power the cycle. In general, the 
sun’s radiation is used to melt a particular PCM at a constant temperature while the 
energy from the subsequent refreezing process is transferred back into the system. Unlike 
sensible heat storage, PCMs enable a more controlled form of energy transfer, as 
discussed by IRENA [8] as the temperature of the medium remains constant throughout 
the phase change. The energy density of PCMs can be much higher than that of sensible 
heat media, especially for metals, according to Ei et al. [17], resulting in reduced 
volumetric requirements. 

Metals as a phase change material are an excellent form of TES. In addition to the 
higher energy density compared to sensible heat media, Kotze and von Backstrom [18] 
discuss how metals typically have a much higher conductivity that greatly enhances the 
heat transfer rate when the TES is in use. The metal selected for the integrated latent heat 
thermal storage in this thesis is Aluminum, based on its thermo-physical properties that 
make it an ideal candidate for Brayton cycle operation. 
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Jiang et al. [19] researched the possibility of using salt as a phase change material, 
exploiting its latent heat abilities, vice sensible heat. As shown in Table 2. the molten 
salt under investigation, Na 2 C 03 -NaCl showed comparable melting temperature and heat 
of fusion specifications to that of standard aluminum. The aluminum, however, has a far 
superior thermal conductivity for heat transfer. 


Table 2. Latent Heat Material Properties of Aluminum and Molten Salt. 

Adapted from [10, 19]. 


Medium 

Melting Temperature 

Heat of Fusion 

Thermal Conductivity 

[C] 

[kj/kg] 

[W/m-K] 

Molten Salt 

637 

332 

2 

Aluminum 

660 

398 

205 
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II. THERMODYNAMIC ANALYSIS 


When designing a new power system, the basic operating limits must first be 
established as a foundation to guide the design of the components. The first step is to 
conduct a thermodynamic analysis of the cycle based on known parameters and 
uncontrolled operating conditions, such as the atmosphere. 

A. GOVERNING EQUATIONS 

The governing equations of mass, energy, and entropy provide a complete 
description the cycle. These equations are coupled with reasonable assumptions of the 
fluid and processes that allow for the calculation of pressure, temperature, heat transfer, 
and work. 


I. Continuity 

The mass of a control volume is accounted for by the continuity equation that 
states that mass can neither be created nor destroyed. The mathematical equation that 
accounts for the sum of all mass entering the control volume, exiting the control volume, 
and any rate of change within, is commonly known. 

dp ^ d{pu) ^ d{pv) ^ g(pw) _Q 
dt dx dy dz 

This equation describes the continuity of mass in Cartesian coordinates. In 
thermodynamics, a broader definition is used to account for the sum of the mass flows 
entering and exiting a system, as well as any change inside. 


dm^ 


dt 


--Z'i’.-Z”- 


( 2 ) 


Power systems such as the cycle analyzed in this thesis, are typically analyzed as 
steady state processes, where the rate of change of mass equals zero. Therefore, the 
continuity equation, as applied to the analysis of a single entry-single exit system, can be 
written as: 
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(3) 


171 = ih = ih 

/ e 

From this expression, the mass flow entering the system equals the mass flow 
exiting the system, and therefore the subscript can be removed. This key relationship will 
further guide the first and second law analysis of the plant. 

2. First Law of Thermodynamics 

The first law of thermodynamics describes the continuity of energy [20]. Similar 
to the law of conservation of mass, energy can neither be created nor destroyed. This law 
accounts for all sources of energy entering and exiting a system: heat, work, and flow. 


dErv ■ ■ . I 1 7 , . 

—— = Qrv ~^rv + nj\ h + —V + gZ -m 

^C.V. C.V. ^ 1 ^ 2 ^ ® ^ I ^ 


h +-V^+gZ 

e 2 ^ ^ e 


(4) 


When analyzing plant cycles at steady state, it is assumed that the change in 
height and the change in velocity are negligible across each process. Additionally, since 
there is no change in the mass flow throughout the cycle, the first law equation is reduced 
to the following format and will be used to describe each process of the cycle. 


(5) 

3. Second Law of Thermodynamics 

The second law of thermodynamics, in conjunction with the first law, provides the 
total basis for a cycle analysis. While the first law conserves energy, the second law 
describes the direction a cycle actually occurs. 


dS^ 


dt 


-ms. -ms +- 


- + S 


( 6 ) 


This equation describes the entropy that exists in the flow as well as heat 
generated in the process. The last term describes the rate that entropy is generated in the 
process. For reversible processes, the generation term is zero but in reality it can never be 
less than zero, thereby showing that entropy is not necessarily conserved. Based on this 
phenomenon, the second law can be reduced to the following form that will be used to 
analyze the power cycle: 
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mis. -s ) + ^^> 0 

B. CYCLE ANALYSIS 

Figure 6 provides a basic illustration of the open air Brayton cycle with the solar 
collector and the connected turbine and compressor. 



Figure 6. Open Air Brayton Cycle with Solar Collector. 

In this cycle, air enters the compressor from the atmosphere (state 1) and is 
compressed to a high pressure gas before entering the solar collector (state 2). The solar 
collector transfers heat to the air, increasing its temperature in a constant pressure process 
(state 3). Finally, the high temperature, high pressure air is used to create work through 
the turbine process before exhausting to the atmosphere (state 4). 

1. Pressure and Temperature 

The initial design of the solar plant is based on the concept of an ideal open 
Brayton cycle. In this cycle, the compressor and turbine processes are assumed to be 
reversible, meaning there is no change in entropy as losses due to friction or other means 
are neglected. The heat transfer occurs by raising the temperature of the fluid on a 
constant pressure line, where state 3 represents an operating limit. 
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A Temperature-Entropy (T-S) diagram of a standard Brayton Cycle, as shown in 
Figure 7, illustrates the change of state of the fluid through each process of the cycle. 
Each number in the diagram represents the corresponding state of the fluid. 


Collector 
■cT 




3 

Turbine 


4 






Comp 


Figure 7. Temperature-Entropy Diagram of Brayton Cycle. Adapted from 

[ 21 ]. 


The working fluid for this cycle is treated as an ideal gas, where the specific heat 
is assumed to be constant. Based on these assumptions, equation (7) can be expanded to 
the following form: 


T P 

s -s.=C ln^-7?ln^ 

e I p rji p 


( 8 ) 


Therefore, the relationship between temperature and pressure for an isentropic 
ideal gas process becomes: 


T P 

C ln^ = y?ln^ 
^ T. P 


(9) 


This expression is simplified by dividing both sides by the specific heat and 
taking the inverse of the equality. 
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(10) 
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_e 

T 

I 


kP.j 


The gas constant for an ideal gas is simply the difference between the specific 
heat for a constant pressure process and the specific heat for a constant volume process 
[ 20 ], 


R = C -C 


P y 


( 11 ) 


The specific heat ratio is the ratio of the specific heats for a constant pressure 
process and constant volume process. 



( 12 ) 


Combining equations (10), (11), and (12), the relationship between temperature 
and pressure for an ideal gas process is defined using the specific heat ratio. 


T 

_e 

T 

i 


( 


V 



r 


(13) 


This relationship is used to define the operating parameters for the compressor 
and turbine processes that are assumed to be isentropic. 

Because the air is exhausted from the turbine to the atmosphere, the pressures at 
the compressor inlet and turbine exit are equal. 


P,=P, (14) 

Assuming the same pressure at these states, and a constant pressure process 
through the solar collector, the pressure ratios of the compressor and turbine are equal. 





(15) 


This relationship between the compressor and turbine is useful because it also 
shows a similar relationship among the temperatures through these processes. 
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(16) 



The known values in this case are the temperature of the air entering the 
compressor (based on atmospheric conditions), and the temperature of the air exiting the 
solar collector, a fixed design point. There are two unknown values: the temperature of 
the air exiting the compressor and the temperature of the air exiting the turbine. Solving 
the equation for the maximum work provides an additional equation that can be used to 
solve for the relationship between these two unknown values. This will produce the 
smallest plant possible to generate the desired amount of power. 

The specific net work for a Brayton cycle is defined as the difference between the 
work of the turbine and compressor. For an ideal gas the difference in enthalpy between 
two states equals the specific heat constant multiplied by the temperature difference. 

'^Net='^T+^C=(K-K) + (K-K) = ^p(T^-T,+T^-T^) ( 17 ) 

The goal of the design is to achieve maximum power density which occurs when 
the derivative of the specific net work is zero. Since the compressor inlet and solar 
collector exit temperatures are fixed constants determined by atmospheric conditions and 
design constraints respectively, the maximum work can be found by varying the 
temperature of the compressor outlet, as shown by Greitzer et al.[22]. 


dw., 


dT^ 



( dT^ 

dT. 

dT 

dp'] 


/ 

dp 1 

c 

3 

-^ + —^ 

2 

= c 

0- 


p 


dT, 

dT, 

~dp] 

p 


dp J 


= 0 


(18) 


After simplification, the relationship between both temperatures is obtained. 



(19) 


From this analysis, the maximum specific net work is achieved when the 
compressor exit temperature and the turbine exit temperature are equal. 


T =T 

4 2 


(20) 
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Combining equation (20) with equation (16), the design temperatures of the 
compressor exit and turbine exit are simply a function of the outside air temperature and 
the specified design temperature of the collector exit. 

( 21 ) 

The specific net work of the cycle can be analyzed as a function of the 
compressor pressure ratio, with the two design points (compressor inlet temperature and 
collector outlet temperature) as scaling factors. By combining equations (13), (17), and 
(20), the work equation becomes: 



In this equation, the pressure ratio includes the exponential factor that is 
determined by the specific heat ratio of the ideal gas. By selecting an ideal gas for the 
cycle and specifying fixed design parameters for the compressor inlet and collector outlet 
temperatures, the work can be plotted against a range of compressor pressure ratios. This 
function will provide a visual representation of the compressor pressure ratio required to 
achieve maximum work. 

The metal selected as the phase change material for the latent heat thermal storage 
is Aluminum. As mentioned before, this is due to its melting temperature and high 
conductivity. As a starting point, the wall temperature of the collector is assumed to be a 
constant temperature, equal to that of the aluminum melting point. The collector air outlet 
temperature was then determined based on a reduced factor of the wall temperature, as 
the air would never actually reach the wall temperature, so long as there is flow. 

= 0.83r = 0.83 (933i<:) = 114.39 K (23) 

The compressor inlet temperature was specified as a fixed value, based on 
average seasonal variation. 


T; =15r = 288.15A' 
17 


(24) 



The cycle uses air as the working fluid that has unique values for the specific heat 
capacity and specific heat ratios, provided in Table 3. 


Table 3. Ideal Gas Properties of Air. 


Property 

Value 

Units 

c 

p 

1.005 

[kJ/kg-K] 

r 

1.41 

[-] 


Using the specified values from equations (23) and (24) while applying the ideal 
gas properties from Table 3. the specific net work of the cycle was plotted as a function 
of the compressor pressure ratio in Figure 8. 



Figure 8. Design Cycle Net Work vs. Compressor Pressure Ratio. 


In this case, where the inlet temperature of the compressor is an assumed fixed 
value, the maximum net work of the cycle corresponds to a compressor pressure ratio of 
5.473. Because the cycle is open, where the turbine exit pressure is atmospheric and there 
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is no assumed pressure drop across the collector, the compressor pressure ratio also 
corresponds to the pressure drop through the turbine. 

The cycle analyzed in Figure 8 is based on an assumed average temperature of the 
air entering the compressor. For the same plant design, the net work was analyzed to 
include more extreme temperatures based on colder and warmer seasons throughout the 
year. Figure 9 illustrates the changes in the cycle net work as a function of the pressure 
ratio. 
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Figure 9. Cycle Net Work vs. Compressor Pressure Ratio for 

Seasonal Variations. 

Based on Figure 9, the cycle achieves a higher maximum net work during times 
when the air temperature is colder that is attributed by the larger temperature differential 
through the solar collector. As the outside air temperature gets colder, the compressor 
pressure ratio required to achieve this maximum work is slightly higher. The computed 



5 10 15 20 25 30 

Compressor Pressure Ratio [N/A] 
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values are provided in Table 4. but in summary the compressor would operate in the 
pressure ratio range of 5-6. 


Table 4. Maximum Specific Net Work and Compressor Pressure Ratios. 


Average Air Temperature Max Work Pressure Ratio 

[C] [K] [kj/kg] [-] 


3 

276.15 

126.446 

5.894 

15 

288.15 

118.520 

5.473 

27 

300.15 

111.007 

5.107 


Equation (13) can be used to solve for the design temperature of the compressor 
outlet, based on the assumed average inlet air temperature. 


T =T 
2 1 




r-1 


1 . 41-1 

288.157^(5.473)^ =A72377K 


(25) 


With the compressor exit temperature known, the temperature at the turbine exit 
can be solved for using the relationship of equation (16). The result shows that the 
maximum work occurs when the compressor outlet temperature equals the turbine exit 
temperature. 


T 



11A39K 


' 288.15j^ ^ 


A12311K = 


(26) 


The lower pressure limit of the cycle is assumed to be atmospheric conditions, 
101.3 kPa. In order to determine the upper pressure limit of the cycle, equation (15) is 
used to solve for the pressure at the compressor exit, knowing that the compressor 
pressure ratio is equal to the pressure drop across the turbine. 


P^=P^=P^^ = l01.3kPa (5.473) = 554.455kPa 
Pi 


(27) 
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2. Specific Work and Heat 

The design specific work of the compressor is determined using the First Law of 
Thermodynamics for a steady state process, adapted for an ideal gas, by using equation 
(5). Compressor work is expressed as a negative value to indicate the work is put into the 
system. 

kJ kJ 

Similarly, the design specific work of the turbine is determined using the 
specified design temperature of the collector exit, and the fact that the turbine exit 
temperature is equal to the compressor exit temperature. Turbine work is expressed as a 
positive value to indicate the work extracted from the system. 

kJ kJ 

w^=h^-h,=c^{T^-T,) = \.m5^—^{llA39K-A12311K) = 202.52'iY (29) 

Finally, the specific heat added to the air via the collector process can also be 
determined using the simple First Law analysis for an ideal gas. For maximum work, this 
value is equal to the specific work of the turbine. 

q^h-h^c (r-r)^ 1.005—11A.59K- 472.377/ir) = 303.523— (30) 

^ ^ /’V 3 2/ ! kg ^ ’ 

The basic design limits and parameters have been determined for the operation of 
the cycle. By establishing an average condition of the air, and setting the design 
temperature limit based on the melting temperature of aluminum, all states were 
determined for the ideal open cycle and provided in Table 5. 


21 



Table 5. Open Brayton Cycle Design Parameters. 


Process 

Temperature 


Pressure 

Specific Work 

Specific Heat 


[K] 


[kPa] 

[kJ/kg] 

[kJ/kg] 



288.150 

P, 

101.325 



Compressor 





-185.148 

- 



472.377 

P. 

554.455 



Collector 


774.390 

P. 

554.455 

- 

303.523 

Turbine 

7 ; 

472.377 

P. 

101.325 

303.523 

- 


3. Efficiency, Power, and Mass Flow 

The thermal efficiency of the system is calculated from the specific net work and 
the specific heat transfer of the collector and corresponds to the maximum work at the 
associated compressor pressure ratio. 


303.523^-185.148^ 

_ kg _^ 

^ ^ 303.523^ 

kg 

Power, defined as the net work rate of the cycle, is a function of mass flow. 

^Ne, = ^^Ne, ( 32 ) 

The initial prototype design is that of a scaled down model, with a design power 
rating of 50 kW. The specific net work, defined as the sum of the compressor and turbine 
work, has been calculated, therefore equation (32) can be used to solve for the design 
mass flow of the cycle. 


( 100 %) = 39% (31) 




m = 


Net 




Net 


50kW 


w 


Net 


Wj + 


kl kJ 

303.523 —-185.148 — 
kg kg 


= 0.422^ 


(33) 


Having determined the design mass flow rate of the air in the cycle, the individual 
plant component powers can now be calculated. 
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kg 


Wc=mwc =0.422^ 
s 


kJ ^ 

-185.148 — 

kg j 


kg 


Wj=mw^ =0A21— 


kJ 

303.523— 

kg 


Q = mq = 0.422^ 
s 


kJ 

303.523— 

kg 


= -l%.2Q4kW 

= 128.204fciy 

128.204fcW^ 


(34) 

(35) 

(36) 


The thermodynamic analysis of the system modeled as an ideal open air Brayton 
cycle, allowed for the calculation of the initial design limits as well as power ratings for 
each of the three components of the system. These calculations provide baseline values 
for sizing the collector as well as selecting a compressor and turbine for the system. Table 
6 . summarizes the resulting design limits. 


Table 6. Summary of Design Limits. 


Parameter 

Value 

Units 

Power 

50.000 

[kW] 

Mass Flow 

0.422 

[kg/s] 

Pressure Ratio 

5.472 

[-] 

ThoI 

774.390 

[K] 
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III. SOLAR COLLECTOR DESIGN 


In Chapter II, a thermodynamic analysis of the cycle yielded a set of design limits 
for the system. These calculations are the basis for designing the first iteration of the solar 
collector’s geometry. 

A. GEOMETRY 

The initial design of the collector is based on the idea of a simple geometry, one 
where reasonable assumptions of heat transfer and flow can be made, as well as the 
ability to support the integration of a phase change metal. 

In this design, the air exiting the compressor enters the solar collector where the 
flow is divided among a stack of evenly sized air ducts. These ducts would be wrapped 
with aluminum phase change metal, so the material is completely integrated into the 
collector heat exchanger. This system takes advantage of the latent heat of fusion. During 
the day, the aluminum is in a liquid state as the sun powers the cycle. When the sun is 
blocked, however, the aluminum undergoes a phase change back to a solid state, but at a 
constant temperature, as it continues to power the cycle. Figure 10 is a side view of the 
cycle that illustrates the TES integration. 


Integrated heat exchanger and thermal storage. 



Figure 10. Side View of Solar Collector with Multiple Stacked Ducts. 
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1 . 


Collector 


Similar to a traditional solar collector receiver, the radiation from the sun is 
directed by an appropriately positioned mirror farm through a small opening or aperture, 
as shown in Figure 11. The radiation is trapped inside the collector with some assumed 
minimal losses to the atmosphere. Heat is then transferred over a larger surface area that 
represents a stack of air ducts. Two design parameters are of particular importance, the 
required area of the aperture and mirror farm. 



a. Aperture 

The radiative heat transfer equation describes the relationship of the area in which 
the radiation passes with the temperatures on each side of the aperture. 

q = -T^) (37) 
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The aperture is a small opening through which the sun’s radiation is captured by 
the mirror farm and then redirected into the collector. The area of the aperture is found by 
rearranging the radiative heat equation. 


‘■Ap 



(38) 


For this case, the collector is modeled as a black body, where the emissivity 
equals 1. The radiation equation is a function of the difference of ambient and internal 
temperatures, both to the fourth power. For the design, the internal temperature is 
modeled as the melting temperature of aluminum, 933 K. As discussed in Chapter II, the 
ambient temperature is assumed to be an average based on seasonal variation, 288.15 K. 
The heat losses through the aperture are assumed to be 5%. Based on these parameters, 
the aperture area is calculated. 


'^Ap 


0.05(1281561T) 


l|^5.67£-8^^j(933"i^" -288.15"i^") 


0.15m" 


In order to provide enough heating to power a 50 kW plant, the area of the 

2 

aperture would be approximately 0.15 m , assuming minimal losses. These losses can be 
customizable. 


b. Mirror Farm 

The area of the mirror farm typically accounts for the bulk of the area required by 
a solar collector. The equation to calculate this area is a function of the reflectivity of the 
mirror that can vary depending on the material of the mirror and particulate collection via 
dust and debris, as well as the losses through the aperture. 




Qap 

pG 


(39) 
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In this case, the reflectivity is assumed to be 85%. The sun’s irradiance is known 

2 

as 1000 W/m . The mirror farm must be large enough to direct enough heat for the cycle 
plus losses from the aperture. 


1.05(128156IT) , 

A... =-k-^ = i58.3lm2 


0.85 


1000 


2 

m y 


For a 50 kW plant, based on reasonable assumptions of the mirror farm and 
considering the losses from the aperture, the area of the mirror farm is approximately 158 


2. Heat Exchanger Ducting 

Having estimated the size of the mirror farm and aperture, the ducting was sized 
to get a complete picture of the collector. In this design, the air temperature is raised to 
the desired design temperature through the straight ducting under uniform heating by the 
by the phase change aluminum. By latent heat phase change, the aluminum maintains 
near melting temperature through the cyclic melting and re-freezing processes. Figure 12 
illustrates the basic model of a single duct. 


Heat 



Figure 12. Single 3-D Square Duct with Molten Aluminum. 
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For the initial design, it was important to estimate the minimum volumetric 
requirements of the aluminum phase change metal required for thermal energy storage. 
To this effect, the boundaries that enclose the aluminum around the air duct were 
considered to be infinitesimally thin rigid walls. In reality, the aluminum would be 
encased by a non-corrosive and preferably conductive material, such as a ceramic with a 
high heat transfer coefficient. The thickness and material properties in this case would 
require additional modeling and would increase the amount of phase change material 
required. 

Since the collector will have a network of stacked ducts, each wrapped with 
aluminum phase change metal, the cross section was chosen based on a geometry that 
supports the reasonable stacking of ducts. A square cross section was chosen for uniform 
heat transfer and flow. 

a. Length 

The heat equation for the internal flow through a volume is a function of the heat 
transfer coefficient, the heat transfer surface area, and the temperature difference of the 
fluid through the duct as it relates to the wall temperature [23]. 

Q = hAAT,^ (40) 

For internal flow convection where the walls are of a constant temperature, the 
heat transfer equation is a function of the log mean temperature difference between the 
inlet and outlet. This relationship accounts for the non-linear temperature rise through a 
duct where the flow temperature approaches the wall temperature at a rate that is 
determined by the mass flow. 


= ■ 


AT - AT 


In 




(41) 


The temperature differences described in this formula are those between the wall 
temperature and the average temperature of the flow at both the inlet an outlet of the duct. 
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( 42 ) 


Ar =r -r, 

o w i 

A7;=r^-r2 (43) 

Substituting equations (42) and (43) into equation (41), the log mean temperature 
difference can then be expanded and simplified to relate the wall temperature with the 
temperatures at the inlet and outlet. 


AT-,,. = 




T -T 


In 


^ 'T' _ T' ^ 


^ 'T' _ T' ^ 


In 




(44) 


The heat transfer coefficient accounts for the flow properties of the air, including 
conductivity. For this model, the flow is treated as laminar and fully developed. For a 
square cross section duct, the Nusselt number is a constant and relates the flow properties 
with the duct cross section. 


Nu^=^ (45) 

k 

The Nusselt number, as shown by Incropera et al. [23], based on the specified 
square cross-section of the duct as well as the assumed flow state, is 3. The heat transfer 
coefficient can now be rewritten as an expression of the Nusselt number, conductivity, 
and hydraulic diameter. 


h = 


Nujjk 


(46) 


The thermal conductivity of the fluid changes as a function of temperature. 
Because the temperature changes non-linearly through the duct, the true conductivity may 
be difficult to model analytically. Therefore, for the initial calculations, a mean 
conductivity of the inlet and outlet was calculated to approximate the thermal behavior 
through the duct. Linear interpolation was used to obtain the true conductivity values at 
the inlet and outlet of the duct. 
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(47) 


k 


kj + ^2 
2 


— kr, 

k = — -^ 


W W 

5 . 609 £- 2 ^^ + 3 . 883 £- 2 ^^ 

m-K m-K 


= 4.746£-2 


m-K 


Since each duct exhibits uniform heating, the total heat transfer area for each duct 
is the sum of the areas of all four sides. 


A=4D,L 


(48) 


The heat equation can then be rewritten by substituting equations (44), (46), and 
(48) in to equation (40). 


Q = 


^ Nuji ^ 
\ J 


(477A) 


T -T 

^3 ^2 


In 

, 1 r -T ,, 

V V w 3 / y 


ANUj,kL[T^-T^) 


In 


_7^ ^ 

\T -T J 
V W 3 / 


(49) 


Using a thermodynamics analysis, the heat transfer can also be described as a 
function of the mass flow, specific heat, and temperature difference through the collector. 


Q = mc^{T,-T,) (50) 

When equations (50) and (49) are combined, the collector length is then solved 
for using the temperature limits identified in Chapter II as well as the flow characteristics 
modeled in this chapter. 


L = 


me, 


ft-r) 


ANu^k (Tj - T^) 


In 


^ 'T' _ T ^ 

^2 

vA. “^3 y 


mc„ 


^In 


ANUjjk 


^ 'T' _ T ^ 

^2 

vA. “^3 y 


(51) 


A summary of the flow parameters and the resulting length for a single duct 
system is provided in Table 7. As mentioned, the wall temperature is modeled based on 
the melting temperature of Aluminum, and the flow was modeled as laminar and fully 
developed, which corresponds to a fixed Nusselt number for the square cross sectional 
duct. 
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Table 7. Collector Duct Heat Transfer Parameters. 


Parameter 

Value 

Units 

Nujj 

2.98 

[-] 

m 

0.422 

[kg/s] 

Kir 

0.047 

[W/m-K] 

K 

1005 

[J/kg-K] 

T 

w 

933.450 

[K] 

T 

-'2 

472.377 

[K] 

T 

774.390 

[K] 

L (Single Duct) 

799.206 

[m] 


Using the values provided in Table 7. the length of a single duct collector was 
calculated at 799.206 m. Since the collector is designed based on a multiple duct stacked 
array, the flow through each duct is proportional to the number of ducts, thus the 
collector length is inversely proportional. The relationship between the collector length 
and the number of ducts in the collector is illustrated in Figure 13. 
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Collector Length [m] 



Figure 13. Collector Length vs. Number of Ducts. 


By stacking ducts, the length of the collector can be dramatically reduced. Some 
examples of collector lengths based on the number of ducts are provided in Table 8. 


Table 8. Collector Length Compared to Number of Ducts. 


Number of Ducts [-] 

Collector Length [m] 

10 

79.921 

50 

15.984 

100 

7.992 

200 

3.996 

500 

1.598 
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b. Hydraulic Diameter 

Sizing the collector length has been useful in analyzing how the overall length can 
be reduced as a factor of N (number of ducts). What also must be known is how the 
length of the collector is related to the hydraulic diameter of the square duct. 

The mass flow equation shows how the flow is related to the cross-sectional area 
of the duct. 


rh = pAu (52) 

For the case of the square duct, the cross-sectional area is the square of the 
hydraulic diameter. 


A = Dl (53) 

When the relationship of equation (53) is substituted back into equation (52), the 
equation still yields two unknowns: the hydraulic diameter and the velocity of the flow. 
The Reynold’s Number for the flow is then used to add an additional equation in order to 
solve for the hydraulic diameter. 


Re = ^ 


(54) 


The Reynold’s Number can then be rearranged as an equation expressing the 
velocity of the flow. 


u = 


Re p 
P^h 


(55) 


Equation (55) can then be substituted in to equation (52) to express the mass flow 
as a function of the Reynold’s Number, viscosity, and hydraulic diameter. 


m = p[Dl) 


yP^kj 


= Re//D^ 


(56) 
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Finally, the hydraulic diameter can then be solved for as a function of mass flow, 
the Reynold’s Number, and the viscosity of the flow. 


D 


h 


m 

Reju 


(57) 


As mentioned, the flow of the air through the duct is modeled as laminar and fully 
developed. For this flow, the Reynold’s Number was assumed to be 2300, the limit 
before the flow transitions to turbulent. 


Similar to the conductivity of the air, the viscosity is also a function of 
temperature that increases non-linearly through the duct for a constant temperature wall 
design. Therefore, a mean viscosity was calculated for the air flow just as it was for the 
conductivity. 

The hydraulic diameter was then calculated based off of the mass flow determined 
from Chapter II as well as the assumptions made for the flow. The flow parameters used 
and resulting size calculations are listed in Table 9. 


Table 9. Collector Duct Flow Parameters and Hydraulic Diameter. 


Parameter 

Value 

Units 

Re 

2300 

[-] 

m 

0.422 

[kg/s] 


3.108E-5 

[kg/s-m] 

A 

5.902 

[m] 

L 

135.413 

[-] 


A single duct collector in this case would have a hydraulic diameter of 5.9 m, but 
just as the collector length, the hydraulic diameter of a single duct also reduces based on 
the number of ducts in the array. For this design, the length to hydraulic diameter ratio 
remains constant at approximately 135, regardless of how many ducts are stacked. 
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B. ALUMINUM REQUIREMENTS 


I. Mass and Volume 

In order to calculate the minimum amount of aluminum phase change metal 
required for the integrated TES, the total amount heat transfer from the collector must be 
determined. 


Q = Qt (58) 

The amount of aluminum required for TES is dependent upon how long the phase 
change metal would be expected to provide latent heat to the plant in the absence of the 
sun. The heat transfer rate from the collector was calculated in Chapter II and is a 
function of the rated power and efficiency of the cycle. 

W 

Q = (59) 

Vth 

Substituting equation (59) for the heat transfer rate in equation (58), the new 
expression for the total heat transfer is expressed as follows: 

W 

Q = (60) 

Vth 

Eor latent heat phase change, the total heat transfer is also a function of the mass 
of the phase change material and its latent heat of fusion. 

In this case, the mass of the Aluminum required is expressed by the density and 
volume. 


’^ai^PaiVai (62) 

Combining equation (60) with (61), and substituting in equation (62), the 
minimum volume of Aluminum is expressed as a function of the rated power, efficiency, 
and the amount of time the phase change metal is expected to provide heating. 


36 



( 63 ) 


W t 

\r — Ner 

^ Al ~ j 

VthP Al^f ,Al 

The Andasol-1 power plant supplies heat from thermal energy storage for up to 
approximately eight hours [12], therefore the design of the aluminum phase change 
volume was determined as a comparable to that standard. Table 10. lists the values used 
to calculate the resulting volumetric requirement using equation (63). 


Table 10. Minimum Aluminum Volumetric Requirement. 


Parameter 

Value 

Units 

w 

Net 

50 

[kW] 

1 

39 

[%] 

Pm 

2712 

[kgW] 

^f,Al 

398000 

[J/kg] 

t 

8 

[hr] 


3.419 

[m'] 


In summary, for a 50 kW rated plant operating within the limits defined in 
Chapter II, the minimum amount of Aluminum required to provide latent heating to the 

■i 

system is 3.4 m . 

2. Collector Duct Thickness 

In the prototype design, each duct within the collector is wrapped in phase change 
aluminum. The minimum amount of Aluminum was determined in the previous section, 
therefore the additional duct thickness from the aluminum can be determined. For 
uniform thickness, the volume can be rewritten in terms of the number of ducts and the 
cross sectional area of the aluminum wrapping, remembering that the length and 
hydraulic diameter are inversely proportional to the number of ducts. 
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For uniformity, the additional thickness can then be determined by subtracting the 
new outer diameter by the hydraulic diameter of the duct, and dividing that quantity by 
two. 



( 66 ) 


For this system, the volume of the aluminum is constant regardless of how many 
ducts are stacked in an array, as it depends on the amount of heating and duration. The 
result, as seen in Table 11. is that the thickness of the aluminum wrapping increases 
when the collector is reduced in length by stacking ducts. 


Table 11. Aluminum Thickness Compared to Number of Ducts. 


Ducts 

[-] 

Collector Length 
[m] 

Hydraulic Diameter 
[mm] 

Aluminum Thickness 
[mm] 

10 

79.92 

590.20 

1.81 

50 

15.98 

118.04 

8.46 

100 

7.99 

59.02 

14.54 

200 

4.00 

29.51 

21.12 

500 

1.60 

11.80 

27.33 
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Table 11. illustrates the disadvantage of spreading a fixed volume of Aluminum 
too thin for a small number of ducts. Stacking as many ducts as possible not only reduces 
the length of the collector, it also provides for a thicker layer of aluminum around the 
ducts, allowing for ease of machining and more uniform heat transfer. 

C. HEAT FLUX 

Of particular importance is the analysis of the heat flux for the prototype collector 
that can be used for sizing optimization. A central solar receiver is designed to focus the 
radiation from the sun on to a relatively small surface area. The expression for the heat 
flux is simply written as the ratio of the heat transfer rate to the heat transfer surface area. 



The goal of central receiver design is to achieve as high of a heat flux as possible. 
Rodriguez-Sanchez et al. [24] showed how current receiver designs can achieve heat 
fluxes on the order of 10 and in some cases 106 W/m . While the prototype solar 
collector proposed is not expected to achieve a heat flux of that magnitude, it is still an 
important factor in size optimization. 

Using MATLAB, the heat flux for the solar collector was first plotted as a 
function of the number of ducts, shown in Figure 14. See Appendix A for the MATLAB 
code. In this case, there was assumed to be no aluminum surrounding the ducts and they 
were stacked as a single array. 
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Heat Flux (W/m ) 



This figure illustrates how the heat flux rises linearly with the number of ducts on 
a logarithmic scale. Based on this scenario the following ducting numbers of Table 12. 
correlate with each approximate heat flux order of magnitude. 


Table 12. Ducting Requirements for Heat Flux Orders of Magnitude. 


Number of Ducts [-] 

Approximate Heat Flux [W/m^] 

50 

10^ 

500 

10^ 

5000 

10^ 


In reality, each duct is wrapped in aluminum. The minimum thickness 


calculations provided in Table 11. showed that the thickness of the aluminum increases 

around each duct as the number of ducts is increased. This leads to a rising total height of 
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the collector, offsetting the decrease in the surface area by reducing the length. Taking 
into consideration the added thickness of the aluminum, the heat flux was then plotted as 
a function of the number of ducts, shown in Figure 15. 



Figure 15. Heat Flux vs. Number of Ducts (With and Without Aluminum). 

Figure 15 provides a more accurate representation of the heat flux behavior based 
on collector sizing, and illustrates how the heat flux actually approaches a limit of 
approximately 2450 W/m . Because of the added thickness around each duct, there is 
eventually a loss of benefit in the collector height sizing that is particularly useful in 
further design refinement. 

D. PRESSURE DROP 

The prototype model presented in this thesis assumed a constant pressure through 

the solar collector, as shown in Chapter II, as a means of conducting an initial cycle 

analysis using a thermodynamic analysis. In reality, there is an expected pressure drop 

across the duct, especially when the length becomes shorter as a result of adding ducts. 
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An analysis of the collector flow was conducted using fluid dynamics as a means of 
quantifying the pressure drop in order to validate the assumption. 

The expression of the pressure drop for a fully developed laminar flow in a non¬ 
circular duct accounts for the wall shear stresses that play a role in the loss of pressure. 




32//Z 


D. 


( 68 ) 


The length to diameter ratio is constant at 135, and the hydraulic diameter in this 
case is a function of N. The maximum velocity is at the exit of the collector where the 
temperature is the highest, and is unknown. Equation (68) can then be rewritten as a 
function of the velocity and N. 


A/^2_3 =2160 


D, 


(69) 


The exit velocity is also found by expanding the mass flow equation in terms of 
an ideal gas. 


P P 

——Au = — —Au^ 
RT^ ' RT^ " 


(70) 


For a constant area duct, the area on each side of the equation will cancel out, as 
well as the gas constant. Assuming a pressure drop across the collector, the outlet 
pressure can be expanded as a subtraction of the pressure drop from the inlet pressure. 
Equation (70) can then be rewritten in terms of the pressure drop. 




^2 y 


(71) 


By combining equations (69) and (71), the velocity at the outlet can be calculated 
by solving the resultant quadratic formula. 


r 

2160 

V 


N 


D 


h J 


/'"s -^ 2^3 +^27^2=0 


(72) 
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The quadratic constants were coded using MATLAB to solve the equation. See 
Appendix A for the MATLAB code used. The inlet velocity was solved as an input array 
based on the scaled mass flow and area as a function of N. 

N 

a = 2160 —u 

b- P 

2 

c = 

2 /T7 2 

^2 

The solution set of the exit velocity as a function of the number of ducts was then 
substituted back into equation (69) and plotted for analysis in Figure 16. 



Figure 16. Collector Pressure Drop vs. Number of Ducts. 

The pressure drop was plotted across a domain of up to 5000 ducts. As the 
number of ducts increased, the pressure drop rose at an exponential rate; however it is 
still very small relative to the collector’s assumed pressure. The collector was assumed to 
operate at a constant 554 kPa, where it would require up to 4000 ducts to see a pressure 
drop of only 1 kPa. Based on this analysis, the assumption of a constant pressure process 
is a valid approximation. 
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IV. OPTIMIZATION 


In Chapter III, the basic sizing parameters of the solar collector for the 50 kW 
prototype plant were established using the operating limits as defined in Chapter 11. The 
power of the plant, however, is a function of mass flow as well as the collector 
temperature rise that is also a function of mass flow. For a fully developed flow in a duct, 
exhibiting a constant wall temperature, the mean temperature of the fluid rises non- 
linearly as it approaches the wall temperature, seen in Figure 17. 



Figure 17. Mean Temperature Profile for Fully Developed Flow, Constant 

Temperature Wall. Source: [25] 

The steepness of this curve depends on the mass flow rate of the fluid. A higher 
mass flow will yield a lower temperature rise, but a lower mass flow yields a lower 
power. This lends to the idea that there is an optimum mass flow where the power is 
maximized as it relates to the temperature rise of the air across the collector. 

An optimization study was conducted analytically to establish an operating curve 
for the plant as a function of mass flow. A similar analysis was conducted using 
numerical methods to capture the real effects of air. These results were compared to see 
how the plant would behave under a range of flows. 
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A. 


ANALYTICAL MODEL 


The temperature of the air exiting the collector was a fixed limit as discussed in 
Chapter II, but in reality, this temperature increases and decreases as a function of mass 
flow. Capturing this phenomenon analytically relies on the relationship between 
thermodynamics and the heat transfer properties of the collector duct. When equations 
(49) and (50) are combined, the logarithmic temperature profile is expressed as a function 
of the flow properties and duct geometry. 


In 




ANUjjkL 

mc^ 


(73) 


The inverse log of both sides of the equation results in an exponential function, 
containing the mass flow term that can be used to scale the temperature profile. 


T -T 
T -T 


exp 


UNu^kL^ 

me „ 


= a = f{m) 


(74) 


For this analysis, the mean conductivity is again assumed to be a fixed average 
although in reality it changes as a function of temperature. The length can be expressed as 
the value calculated in Chapter III for a single duct collector, as the number of ducts that 
scale both the length and flow offset each other. The collector inlet temperature is then 
expressed in terms of the scaling function, wall temperature, and outlet temperature. 


T,={l-a)T^ + aT, (75) 

As discussed in Chapter II, from equation (21) the maximum work occurs when 
the collector inlet temperature equals the square of the compressor inlet temperature and 
collector outlet temperature. 


Substituting equation (21) into equation (75) produces a new expression for the 
collector outlet temperature that is assumed to be unknown. 
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(76) 


T^={\-a)^ ^ + 2[a-a^)^^ + a 


T r, 2 t; 

W 3 , 3 


By collecting like terms, equation (76) can then be rearranged in the form of a 
second order polynomial equation, where the scaling function contains the mass flow 
term. 


r/ + 


2 


T 




w 



(l-«) 2 

a 


(77) 


Equation (77) is then solved for the collector outlet temperature for each value of 
mass flow across a specified domain using the quadratic formula and associated quadratic 
constants. 


T = 


-b±\]b^ -4ac 


2a 


(78) 


a = l 


b = 





2 

-1 

T - 


7 

a 


The quadratic formula was solved for the collector outlet temperature as a 
function of mass flow. See Appendix B for the code executed in MATLAB to solve the 
equation. 

Having solved for the collector outlet temperature, the plant power can then be 
calculated using equation (32) for ideal gases that expresses the net work in terms of 
temperature states. 


Ke, =H^T+Wc) = ^C^[{T^-T,)+{T^-T^)] 

Recalling again that for maximum work, the compressor outlet temperature and 
the turbine outlet temperature are equal for the design prototype of the smallest size, the 
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power equation can be rewritten in terms of the compressor inlet temperature and 
collector outlet temperature that was just solved. 




(r3-7^) + (7j-7^) ='”c^(7’3+^l-2^/^) (79) 


The compressor inlet temperature again is a fixed uncontrollable constant based 
on the atmospheric conditions, but similar to the methods of Chapter II, it can be assumed 
based on an average of seasonal variation. Substituting the fixed limits into equation (79), 
as well as the solution for the collector outlet temperature, the power can then be plotted 
as a function of mass flow. The resulting operating curve for up to 2 kg/s of flow is 
shown in Figure 18. 



Figure 18. Analytical Solution of Power vs. Mass Flow for a 50 kW Plant. 


The analytical solution of the plant power as a function of mass flow yields a 
maximum power of the 50 kW plant design. Based on the operating curve, the previously 
calculated mass flow from Chapter II correlates with the design power of 50 kW on the 
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curve. The analysis, however, shows that the power can be maximized by inereasing the 
mass flow to nearly 1 kg/s. Table 13. provides a summary of the exaet calculations for 
both design and maximum power eonditions. 


Table 13. Design and Maximum Power Calculations for a 50 kW Plant. 

Mass Flow [kg/s] Power [kW] 

Design 0.42 50.00 

Maximum 0.96 58.18 


The analytieal analysis represents the potential to exceed design power by 16%. 
In this case, the lowered eollector temperature differenee is offset by the increased flow 
through the duet. Figure 19 shows how the eollector temperature difference deereases in 
an exponential decay as the mass flow in increased. 



Figure 19. Analytical Solution of Collector Temperature Difference vs. Mass 

Flow for a 50 kW Plant. 
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The initial design prototype produced a basic operating collector temperature 
difference of approximately 301 K, but that difference can drop as low as 192 K when 
maximum power is preferred. The analysis shows that the plant could support a 
significant increase in flow relative to design for more power. 

The plant efficiency was also plotted against the same domain of mass flow to 
observe the difference between design and maximum power efficiencies, as shown in 
Figure 20. 



Figure 20. Analytical Solution of Thermal Efficiency vs. Mass Flow for a 50 

kW Plant. 

From Chapter II, the thermal efficiency of the plant modeled at 50 kW was 
calculated at 39%, but at maximum power, the thermal efficiency drops to 31%. This 
observation is particularly useful for the operational considerations of the plant. Solar 
thermal plants aim to maximize power over efficiency, so a higher mass flow would be 
preferred, especially during the day when the plant is powered by the sun. However, 
during night time operations or inclement weather conditions, the thermal energy stored 
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by the phase change aluminum must be conserved. The reduced flow in this case is 
preferred to minimize the power in a way that meets demands up to design limits but 
maximizes the duration of the phase change thermal storage capacity. 

Finally, the heat exchanger effectiveness was analyzed to observe the 
performance of the collector across the operating range of mass flows. The heat 
exchanger effectiveness is expressed as a percentage ratio based on the maximum heat of 
the collector. 


E - 


q 


q 


Max 


(80) 


The maximum heat from the collector is achieved when the collector outlet 
temperature equals the wall temperature. For the constant wall temperature duct model, 
the effectiveness of the heat exchanger is maximized when the mass flow is nearly zero, 
which is not ideal as the power is then minimized. Similar to the efficiency analysis, the 
heat exchanger effectiveness drops as a function of mass flow, as illustrated in Figure 21. 



Figure 21. Analytical Solution of Heat Exchanger Effectiveness vs. Mass Flow 

for a 50 kW Plant. 
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For the design parameters, the heat exchanger effectiveness is 0.65, but drops to 
0.37 when maximum power is achieved. This operating range shows how the heat 
exchanger should perform during day or night operation to maximize the use of the 
molten aluminum thermal storage. 

B. NUMERICAL MODEL 

The analytical model was useful in determining approximate curves to describe 
the behavior of the plant and solar collector. As stated previously, some broad 
assumptions were required to simplify the optimization model in a way that could be 
solved analytically. The flow was modeled as laminar and fully developed. The 
corresponding Nusselt number for a cross-sectional duct could then be used to define the 
heat transfer coefficient in the heat equation for the constant temperature wall duct. 
Additionally, both the thermal conductivity and viscosity of the flow were assumed to be 
fixed constants, calculated as mean values based on the inlet and outlet temperatures of 
the duct. In reality, as previously discussed, these values vary with temperature, which 
can change the profile through the duct. 

In order to better understand how the real flow through the duct would behave, a 
model was constructed to simulate similar conditions to the design while removing the 
assumptions of the temperature dependent flow parameters. ANSYS-CFX software was 
used to numerically solve for the flow behavior, using a model constructed with 
SolidWorks. 

1. Model Construction 

A single duct was created using SolidWorks software. By modeling a single duct, 
the behavior can be scaled accordingly for a multiple duct system. The length to 
hydraulic diameter ratio remains constant regardless of how many ducts are stacked in 
the solar collector. Therefore, the dimensions of the duct are irrelevant, so long as the 
geometric ratio is 135, and the mass flow rate through the duct is scaled accordingly. 
Figure 22 provides an illustration of the model, which represents the air duct under 
consideration. 
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Figure 22. SolidWorks Drawing of Single Duct Model. 



For the prototype duct, the length and hydraulic diameter were specified at 400 m 
and 2.95 m respectively. This model would correspond to a two-duct collector of a 50 
kW plant. The geometry was then imported into ANSYS-CFX software for a full-scale 
analysis of the flow. 

The first step of the analysis was to create a working mesh by which the model 
can be numerically solved. The resulting generated mesh, as shown in Figure 23 has a 
total of 28764 elements and 39200 nodes (see Appendix C). 



Figure 23. ANSYS-CFX Mesh for Single Duct Model. 
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2 . 


Boundary Conditions 


The next step of the solution process was to establish the fluid domain and 
boundary conditions for the model. The fluid model for this analysis was selected such 
that the air behaved as an ideal gas and that the flow was laminar, but incorporated real 
effects in order to achieve more realistic results. In doing so, the effects of changing 
viscosity and conductivity were included. This real effect of the flow could not be easily 
captured analytically, but is accounted for numerically using Sutherland’s law that is 
expressed mathematically for both flow properties. 


( 


3 


(P = (Pref 


T 

\refj 


,+S 

ref _ 

T + S 


(81) 


In this formula, the flow parameter changes as a function of temperature relative 
to a reference point. The “S” refers to Sutherland’s constant and is represented in units of 
Kelvin. By selecting this option, the program ANSYS-CFX invokes the relationship 
throughout the solution process. 


Table 14. Sutherland’s Constant and Reference Properties. 


Property 

Value 

Units 

5 

110.4 

[K] 

T 

^ref 

273.15 

[K] 

^ref 

1.72(10'^) 

[kg/m-s] 

^ref 

2.43(10'^) 

[W/m-K] 


The duct has a total of six boundaries: inlet, outlet, and four walls. The boundary 
conditions at the inlet and outlet were specified as mass flow, scaled for a two-duct 50 
kW model (0.211 kg/s). In order to resolve continuity, however, a temperature or 
pressure must be specified at the inlet mass flow boundary. The inlet temperature in this 
case was first specified as the same temperature calculated in Chapter II (472.377 K). 
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Each wall boundary for this model was established as a constant temperature wall, 
equal to the melting temperature of Aluminum (933 K). These wall boundaries were also 
expressed as no slip boundaries, where the velocity of the flow at the wall is zero. Figure 
24 illustrates the default domain and boundaries for the fluid model. 



1000C 30 000 


Figure 24. ANSYS-CFX Default Domain with Boundaries for Single Duct 

Model. 

The model had a convergence target of 10-^ with the maximum number of 
iterations set to 500 with automatic time stepping. 

3. Initial Run 

The model was run once to check the flow behavior and properties at the outlet 
compared to the analytical calculations. Figure 25 is the plot of the convergence of mass 
and momentum for the initial run. 
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Figure 25. ANSYS-CFX Mass and Momentum Solution Convergence for 

Single Duct Model. 

According to the plot, the mass and momentum equations were resolved with 
good convergence. The results of the simulation were then observed to check the flow 
properties at the duct outlet. 

The first parameter checked was the temperature of the air at the outlet of the duct 
to compare with the value calculated in Chapter II. Figure 26 illustrates the temperature 
profile at the outlet of the duct. 
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Temperature 

Temp 

m 9.330e+002 


8.157e+002 


6.984e+002 


5.811e+002 

^ 4.638e+002 
[K] 



Figure 26. Collector Duct Outlet Temperature Profile. 


The average temperature at this location was calculated to be 742 K. This value 
was slightly smaller than the temperature upper limit defined in Chapter II, which was an 
approximation based on an assumed reduction from the wall temperature (772 K). The 
numerical solution indicates a 20% reduction as opposed to the assumed 17% reduction 
previously defined, which shows that the analytical assumptions are reasonable. 

An important assumption about the flow in the analytical modeling process was 
that of the Reynold’s Number that was fixed at 2300, representing the limit before a flow 
transitions to turbulent. In order to validate that assumption, the Reynold’s Number was 
checked at the outlet of the duct where the velocity was expected to be a maximum. 
Figure 27 shows the velocity profile of the flow at the duct outlet. 
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Figure 27. Collector Duct Outlet Velocity Profile. 


Z 



The Reynold’s Number was computed using an expression in CFX, based on the 
average velocity of the flow at the outlet and equals 2620 for this model. This value is 
slightly higher than the assumed value, and just inside the transition region, but again 
shows that the analytical assumptions were reasonable. 

In Chapter III, the pressure drop in the collector was analyzed analytically as a 
function of the number of ducts. Based on that analysis, the pressure drop increased 
exponentially as a function of N. For collectors with two ducts, the pressure drop was 
found to be 2.46(10''^) Pa (see Appendix A). The pressure drop in the numerical model 
was quantified by observing the minimum and maximum pressures of the flow in the duct 
that correlate with the inlet and outlet of the duct, respectively, and was found to be 
4.60(10’'^) Pa (see Appendix C). 

4. Optimization 

Having established satisfactory results from the initial simulation run, the model 
was then analyzed using the Optimization tool in ANSYS, (see Appendix D). This 
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feature allows the user to specify a range of design inputs, in this case mass flow, instead 
of just one. The software runs multiple simulations across evenly spaced values in that 
range and stores the resulting output parameter(s) specified by the user. For this model, 
the three output parameters of significance were the collector outlet temperature, the 
collector temperature difference, and the resulting plant power. 

a. Collector Outlet Temperature 

In order to capture a similar comparison numerically, the inlet temperature of the 
collector duct boundary was specified as an expression, using the same relationship 
defined by equation (21) for maximum work. This resulted in ANSYS-CFX calculating 
the inlet temperature iteratively throughout the solution process. After 25 simulations, the 
outlet temperature was recorded (see Appendix E) and plotted in excel with the analytical 
solution. Figure 28 illustrates the comparison of the outlet temperature as a function of 
mass flow between the analytical and numerical models. 


1000 



400 

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 

Mass Flow [kg/s] 

Figure 28. Comparison of Analytical and Numerical Solutions for Collector 

Outlet Temperature vs. Mass Flow. 
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Based on the plot, the results between both analysis methods were very similar, 
with the maximum error occurring at a mass flow value of approximately 0.1 kg/s. 

b. Collector Temperature Difference 

Similarly, the optimization tool was used to establish a curve for the collector 
temperature difference as a function of mass flow. Figure 29 illustrates comparison of 
both analytical and numerical models. 

450 
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Figure 29. Comparison of Analytical and Numerical Solutions for Collector 

Temperature Rise vs. Mass Flow. 

The results were again very similar in this case, with the maximum error also 
occurring at a mass flow value of approximately 0.1 kg/s. 

c. Power 

Most importantly, the plant power was observed using the optimization tool to 
compare with the analytical model when more real effects of the flow were invoked. In 
this case, the output parameter, power, was expressed as equation (79), such that the 
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numerical solution of the power only relied on the compressor inlet temperature (fixed) 
and the numerical solution of the collector outlet temperature as a function of mass flow. 
Figure 30 shows the resulting comparison between both the analytical and numerical 
curves. 


Analytical 



0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 

Mass Flow [kg/s] 

Figure 30. Comparison of Analytical and Numerical Solutions for Power vs. 

Mass Flow. 


The analytical model, which assumed a constant average viscosity and 
conductivity, yielded a maximum power of 58 kW corresponding to a mass flow of 
approximately 1 kg/s. For the numerical model, the optimization showed that the 
maximum power of 59 kW could actually be achieved at a higher mass flow, 
approximately 1.26 kg/s. In general, the curve shapes were very similar, with the 
calculated maximum power values having only an approximate 1 kW difference. Table 
15. shows a summary of the exact calculations from each model. 
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Table 15. Maximum Power and Pressure drop Based on Analytical and 

Numerical Solutions 


Model 

Mass Flow 

Maximum Power 

Pressure Drop 

[kg/s] 

[kW] 

[lO"* Pa] 

Analytical 

0.96 

58.18 

2.46 

Numerical 

1.26 

59.46 

4.62 


In summary, the numerical model was intended to provide a more exact solution 
to compare against the analytical solution, which made more broad assumptions of the 
flow. The comparison of both models indicated that the assumptions and approximations 
made to enable an analytical solution were fairly accurate and could be used for 
extremely fast approximations. 
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V. DISCUSSION 


Thermal Energy Storage systems in use today are very complex and require the 
pumping and circulation of multiple fluids to power solar thermal plants during non¬ 
daylight hours. These complexities can add additional operational burden, increased plant 
maintenance, and drive down overall plant efficiency. The purpose of this thesis was to 
establish a new design centered on simplicity, completely removing any pumping of hot 
liquid to foster the automatic transition between solar power supply and thermal energy 
storage. Many plants today operate on a Rankine cycle, but the prototype design in this 
thesis is based off of a Brayton cycle, so that air is the only moving fluid, further 
simplifying the plant. 

In order to only have only one moving fluid for this plant, the thermal energy 
storage integrated into the collector for this design utilized the latent heat of phase change 
of aluminum as this minimized the volumetric requirements of the material. Other metals 
may be possible but Aluminum was suitable for operation when a Brayton cycle is used, 
because of its ready availability and relative low cost. 

The prototype design model presented in this thesis intended to establish some 
basic operating limits and sizing specifications, as well as minimum volumetric 
requirements of the aluminum phase change metal. The assumptions of the collector flow 
were the largest contributor to the geometric design. For the 50 kW prototype, the 
Reynold’s Number was assumed to be 2300 based on those assumptions. As stated in 
Chapter III, the hydraulic diameter is a function of mass flow, which scales 
proportionally with plant power. If the prototype were scaled up, the resulting collector 
height would not be feasible. A 50 MW plant, for example, would have a minimum 
height of 5900 m. In order to reduce the collector height for a larger plant, the flow 
should not be modeled as laminar. By increasing the Reynold’s Number to turbulent 
flow, the minimum collector height as well as the length could be dramatically reduced. 
Another option would be to change the material to a higher temperature metal. 
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The model intended to estimate the minimum volumetric requirements of 
aluminum to sustain up to eight hours of power without the sun. This calculation required 
that the boundaries encasing the aluminum were infinitesimally thin and rigid. In reality, 
the aluminum would be wrapped with a non-corrosive material, preferably of a high 
conductivity to maximize the heat transfer. 

Solar central receivers in operation and under design achieve extremely high heat 

fluxes. The prototype design in this thesis achieved a maximum heat flux on the order of 
-1 

10 , limited by the aluminum wrapping. The analysis in Chapter III showed how the 
thickness of the aluminum around the ducts increases as the ducts are stacked, although 
the length and hydraulic diameters are reduced by a factor of N. The limit assumed some 
customizable losses from the aperture and mirror farm, which can be explored further to 
seek ways to improve the heat flux for this design. 
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VI. CONCLUSION 


The design is a new concept; therefore, a basic thermodynamic analysis was used 
to conduct a complete cycle analysis of the system to establish operating limits necessary 
to design the collector. The upper temperature limit was established based on the melting 
temperature of Aluminum. By fully integrating latent heat thermal storage into the 
collector, the efficiency loss would be very minimal, limited only to heat losses to the 
atmosphere. These losses are easily mitigated with proper insulation technologies. 

The collector was sized based on the limits defined from the cycle analysis. The 
collector duct geometry was selected based on the idea of stacking multiple ducts, which 
reduces the length and hydraulic diameters as the ducts are stacked in a column. Stacking 
ducts also increases the heat flux by minimizing the heat transfer surface area in the 
collector. The model was assumed to have no pressure drop, which was validated by 
analytical and numerical methods, that showed very little pressure drop. 

The plant behavior as it pertains to mass flow was analyzed using both analytical 
and numerical methods. Using the equations from thermodynamics and heat transfer, the 
plant power and performance characteristics were calculated and plotted as a function of 
mass flow to see how these parameters behave under different conditions. The results 
showed that the design power could actually be exceeded by increasing the mass flow to 
a certain extent. Other performance characteristics were also analyzed to establish 
operating bands of the plant between day time and night time. The power and 
temperature calculations were compared to a numerical model which incorporated more 
real effects of the air and flow. The results from both models showed some similarities, 
and the assumptions made in the analytical model were validated as good 
approximations. 
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VII. RECOMMENDATIONS 


The research conducted in this thesis has hopefully established a foundation on 
which the prototype plant can be further developed. In an effort to guide the continued 
effort in this design study, the following recommendations are provided: 

1. Analyze larger scale plant design 

2. Revise flow model for larger collector sizing 

3. Investigate suitable materials to encase the phase change metal around the 
ducts, such as high conductivity ceramics 

4. Investigate modifications to the cycle such as a bottoming Rankine cycle 
to improve the plant efficiency 
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APPENDIX A. MATLAB CODE EOR ANALYSIS OE PRESSURE 
DROP AND HEAT ELUX OE 50 KW PLANT 


% LT Collin Roof, Thesis 

%% Collector Pressure Drop and Heat Flux 


Power = 50000; 
nu = 0.39; 

Q = Power/nu; 

PI = 101325; 

T1 = 288.15; 

Tw = 933; 
f = 0.83; 

T3 = f*Tw; 

T2 = sqrt(Tl*T3 
gamma = 1.4; 

P2 = PI*(T2/T1) 


% [W] Rated Power 
% [N/A] Thermal Efficiency 
% [W] Collector heat xfer 
% [Pa] Compressor inlet pressure 
% [K] Compressor inlet temperature 
% [K] Collector wall temperature 
% [N/A] Collector outlet reduction 
% [K] Collector outlet temperature 
; % [K] Collector inlet temperature 

% [N/A] Air specific heat ratio 
(gamma/(gamma-1)); % [Pa] Collector inlet pressure 


mu3 = ( (1.458e-6)*T3^(3/2))/(T3 + 110.4) ; 

% (Sutherland's Law) 

R = 287; % [kJ/kgK] Air Gas Constant 

D = 5.902; % [m] Single Duct Diameter 

L = 799.206; % [m] Single Duct Length 

C = L/D; % [N/A] 


[kg/ms] Collector out viscosity 


N = 5:50:5000; 
mdot = 0.422./N; 
d = D./N; 

1 = L./N; 

A = d.^2; 

As = N.*1.*d; 

V2 = (R*T2/P2).*mdot./A; 

Qf = Q./As; % [W/m2] 


% [N/A] Number of Ducts 
% [kg/s] Mass Flow/duct 
[m] Diameter/duct 
[m] Length/duct 
[m^2] Area/duct 

% [m^2] Heat xfer surface area 

% [m/s] Inlet velocity/duct 
Heat flux 


% Quadratic Constants 
a = 16*C*(mu3./d); 
b = -P2; 

C = P2*(T3/T2).*V2; 


% Solve Quadratic Equation for V3 (Collector Outlet) 
for i = 1:length(N) 

V3(:,i) = roots([a(i) b c(i)]); 
end 


% Removing erroneous root 
V3 = V3(2, : ) ; 

% Plug V3 back into formula for Pressure Drop 

dP = ((16*C*mu3)./d).*V3; % [Pa] Collector Pressure Drop 

% Plot Pressure Drop vs. Number of Ducts 
figure(01) 

set(gca, 'fontsize' ,15, 'linewidth', 2) 
hold on 

plot(N, dP, 'LineWidth' ,4) 
xlabel( 'Number of Ducts [-]' ) 
ylabel (' Collector Pressure Drop [Pa]') 
grid on 

% Plot Heat Flux vs. Number of Ducts 
figure(02) 

loglog(N, Qf, 'Linewidth' ,4) 
set (gca, 'fontsize' ,15) 
hold on 

xlabel( 'Number of Ducts [-]' ) 
ylabel ('Heat Flux [W/m^2]' ) 
grid on 
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APPENDIX B. MATLAB CODE EOR ANALYTICAL 
OPTIMIZATION OE 50 KW PLANT 


%% LT Collin Roof, Thesis, Solar Plant Optimization, 50kW Plant 
%% Introduction 

% The following script optimizes mass flow of the solar collector Brayton 
% cycle for maximum power. 

% Assumptions: 

% 1) Air - Ideal Gas 

% 2) Duct(s) - Constant area square cross section, L/a = constant 
% 3) Fully developed flow, Nu provided from text. Re = 2300 
% 4) Average constant conductivity through duct 
% 6) Ideal, Open Cycle - Compressor and Turbine processes 
% reversible/adiabatic 
% 7) Design Power - 50 [kW] 

clear all 
close all 
clc 

%% Setup 

T1 = 288.15; % [K] Compressor Inlet Temperature 

PI = 101.325; % [kPa] Compressor Inlet Pressure 
cp = 1.005; % [kJ/kg-K] Specific Heat Capacity of Air 

gamma = 1.41; % [N/A] Specific Heat Ratio of Air 

R = 0.287; % [kJ/kg-K] Gas Constant of Air 

Tw = 933.45; % [K] Wall Temperature (Melting temperature of A1... 

% and max operating temperature of cycle) 

P4 = PI; % [kPa] Turbine Outlet Pressure 

Nu =2.98; % [N/A] Nusselt Number for Square X-Section Duct 

k = 4.746e-5; % [kW/m-K] Average Conductivity of Air through duct 
L = 799.206; % Length of a single duct 

%% Initialization 

% Set domain for Mass Flow 

mdot = linspace(0.01,2,135); % [kg/s] Mass Flow 

% Coefficient function of mdot from heat equation 
alpha = exp((4*Nu*k*L/cp)./mdot); 

% Quadratic Solver for T3 using max power relationship with T2 
a = 1 ; 

b = ( (2.*(alpha-alpha.^2).*Tw./T1)-1)./((alpha.^2)./T1) ; 
c = (((1-alpha) .^2) .*(Tw^2) ./Tl) ./( (alpha.^2) ./Tl) ; 

% Root solver loop 
for i = 1:length(mdot) 

T3(:,i) = roots ( [a b(i) c(i)]); 
end 

% Removing 2nd root array as not realistic 
T3 = T3 (1, : ) ; 


% Solving for T2 using max power relationship 
T2o = sqrt(T1.*T3); % Optimizer T2 (Varies with T3) 

% Delta T through duct 
dT = T3-T20; 

% Cycle Power as a function of Mass Flow, Max Power 
P = mdot.*cp.*T1.*((T3./Tl)-2.*sqrt(T3./Tl)+1); 

[MI] = max (P) ; 
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% Heat Exchanger Effectiveness 
epsilon = (dT)./(Tw - T2o); 

% Cycle Efficiency 
eta = P./(mdot.*cp.*dT); 

%% Results/Figures 

% Delta T Through Duct as a Function of Mass Flow 
figure(01) 

set (gca, 'fontsize' ,15, 'linewidth', 2) 
hold on 

plot(mdot, dT(l,:), 'LineWidth' ,4) 

plot(0.422, dT(29), '.g', 'markersize' , 30) 

hold on 

plot(mdot(I), dT(I), '.r', 'markersize', 30) 
hold on 

text(0.422, dT(29), sprintf([' Design' ]),.. . 

'VerticalAlignment' , 'bottom' , 'HorizontalAlignment' , 'Left' , ... 

'fontsize' ,15) 
hold on 

text(mdot(I), dT(I), sprintf([ 'Maximum' ]),.. . 

'VerticalAlignment' , 'bottom' , 'HorizontalAlignment' , 'Left' , ... 

'fontsize' ,15) 
xlabel( 'Mass Flow [kg/s]' ) 
ylabel (' Delta T [K]' ) 

%title('Collector Delta T vs. Mass Flow') 
grid on 

% Power as a function of mass flow 
figure(02) 

set (gca, 'fontsize ', 15, 'linewidth', 2) 
hold on 

plot(mdot, P(l,:), 'Linewidth' ,4) 

hold on 

plot(mdot(I), M, '.r', 'markersize', 30) 
hold on 

text(mdot(I), M, sprintf ([' \nMaximum:\nPower = 58 [kW]\nMass Flow = 1 [kg/ 
s] ']),.. . 

'VerticalAlignment' , 'top' , 'fontsize' ,15) 

hold on 

plot(0.422, 50, '.g', 'markersize', 30) 

text(0.422, 50, ... 

sprintf (' \nDesign:\nPower = 50 [kW]\nMass Flow = 0.4 [kg/s]'),... 

'VerticalAlignment' , 'top' , 'fontsize' ,15) 

hold on 

xlabel( 'Mass Flow [kg/s]' ) 
ylabel( 'Power [kW]' ) 

%title('Cycle Power vs. Mass Flow') 
grid on 

% Heat exchanger effectiveness as a function of mass flow 
figure(03) 

set(gca, 'fontsize ', 15, 'linewidth', 2) 
hold on 

plot(mdot, epsilon, 'Linewidth' ,4) 
hold on 

plot(0.422, epsilon(29), '.g', 'markersize', 30) 
hold on 

plot(mdot(I), epsilon{I), '.r', 'markersize', 30) 
hold on 

text(0.422, epsilon(29), sprintf([ 'Design' ]),.. . 

'VerticalAlignment' , 'bottom' , 'HorizontalAlignment' , 'Left' , ... 

'fontsize' ,15) 
hold on 

text(mdot(I), epsilon(I), sprintf([ 'Maximum' ]),.. . 

'VerticalAlignment' , 'bottom' , 'HorizontalAlignment' , 'Left' , ... 

'fontsize' ,15) 
xlabel( 'Mass Flow [kg/s]') 
ylabel (' \epsilon [-]' ) 
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%title('Heat Exchanger Effectiveness vs. Mass Elow') 
grid on 

% Cycle efficiency as a function of mass flow 
figure(04) 

set(gca, 'fontsize' ,15, 'linewidth', 2) 
hold on 

plot(mdot, eta, 'LineWidth' ,4) 
hold on 

plot(0.422, eta(29), '.g', 'markersize' , 30) 

hold on 

plot(mdot(I), eta(I), '.r', 'markersize', 30) 
hold on 

text(0.422, eta(29), sprintf([ 'Design' ]),.. . 

'VerticalAlignment' , 'bottom' , 'HorizontalAlignment' , 'Left' , ... 

'fontsize' ,15) 
hold on 

text(mdot(I), eta(I), sprintf([' Maximum' ]),.. . 

'VerticalAlignment' , 'bottom' , 'HorizontalAlignment' , 'Left' , ... 

'fontsize' ,15) 
xlabel( 'Mass Flow [kg/s]' ) 
ylabel( '\eta [-]' ) 

%title('Cycle Efficiency vs. Mass Flow') 
grid on 

%% Temperature Distribution Through Single Duct, Fixed Mass Flows 

% Analyzing the Temperature Profile Through the Duct for Mass FLows 

T2 = 472.491; % [K] Collector Inlet Temperature 

mdotl = [0.01 0.42 0.96 1.5 2]'; % [kg/s] Mass Flow Samples 

X = linspace(0,800,1000); % [m] Duct Length (Dependent Variable) 

1 = x/800; % [N/A] Normalized Length 

alpha = exp(((4*Nu*k/cp)./mdotl).*x); % [N/A] Scaling Function 
T = (1-1./alpha)*Tw + (1./alpha)*T2; % [K] Temperature Across Duct 
% Results/Plot 
figure(05) 

set (gca, 'fontsize' ,15, 'linewidth', 2) 
hold on 

plot(l,T(l,:), 1,T(2,:), '-g' , 1,T(3,:), '-r', 1,T(4,:), 1,T(5,:),... 

'Linewidth' ,4) 

xlabel( 'Normalized Duct Length [-]' ) 
ylabel( 'Average Temperature [K]' ) 

h=legend ('0.010 [kg/s] ' , ' 0.422 [kg/s] (Design)' , . .. 

'0.960 [kg/s] (Maximum) ' , '0.150 [kg/s] ' , . . . 

'0.200 [kg/s]' ) 
set(h, 'Position',... 

[0.219433992346939 0.572591553287982 0.286607142857143 0.204761904761905]) 
%title('Average Temperature Distribution Through Single Duct') 
grid on 
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APPENDIX C. ANSYS-CFX SOLUTION REPORT FOR 50KW 

PLANT 


/\NSYS 


Title 

Solar Collector Single Duct 

Author 

LT Collin Roof 

Date 

2017/05/23 20:45:37 


Contents 

1. File Report 

Table 1 File Information for CFX 1 

2. Mesh Report 

Table 2 Mesh Information for CFX I 
Table 3 Mesh Statistics for CFX 1 

3. Physics Report 

Table 4 Domain Physics for CFX 1 
Table 5 Boundary Physics for CFX 1 

4. Solution Report 

Table 6 Boundary Flows for CFX 1 
Table 7 Forces and Torques for CFX 1 

5. Pictures 

Figure 1 Pressure Profile 
Figure 2 Temperature Profile 
Figure 3 Velocity Profile 
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1. File Report 

Table 1. File Information for CFX 1 


Case 

CFX 1 

File Path 

H:\MyDocs\Thesis\Model\SingeDuct(Scaled)_files\dpO\CFX-l\CFX\Fluid Flow 

CFX_061.res 

File Date 

22 May 2017 

File Time 

10:25:13 PM 

File Type 

CFX5 

File 

Version 

17.2 


. Mesh Report 

Table 2. Mesh Information for CFX 1 


Domain 

Nodes 

Elements 

Tetrahedra 

Wedges 

Pyramids 

Hexahedra 

Polyhedra 

Default Domain 

39200 

28764 

0 

0 

0 

28764 

0 


Table 3. Mesh Statistics for CFX 1 


Domain 

Minimum 
Face Angle 

Maximum 
Face Angle 

Maximum Edge 
Length Ratio 

Maximum 
Element Volume 
Ratio 

Connectivity 

Range 

Default 

Domain 

90 [ degree ] 

90 [ degree ] 

1.01691 

1.00006 

1 

8 


76 


























3. Physics Report 

Table 4. Domain Physics for CFX 1 


Domain - Default Domain 

Type 

Fluid 

Location 

B30 

Materials 

Air Ideal Gas 

Fluid Definition 

Material Library 

Morphology 

Continuous Fluid 

Settings 

Buoyancy Model 

Non Buoyant 

Domain Motion 

Stationary 

Reference Pressure 

5.5491e+02 [kPa] 

Heat Transfer Model 

Total Energy 

Include Viscous Work Term 

Off 

Turbulence Model 

Laminar 


Table 5. Boundary Physics for CFX 1 


Domain 

Boundaries 


Boundary - Inlet 


Type 

INLET 


Location 

Inlet 


Settings 


Flow Direction 

Normal to Boundary Condition 

Default Domain 

Flow Regime 

Subsonic 


Heat Transfer 

Static Temperature 


Static Temperature 

T2 


Mass And Momentum 

Mass Elow Rate 


Mass Flow Rate 

MassElow 


Mass Flow Rate Area 

As Specified 
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Boundary - Outlet 

Type 

OUTLET 

Location 

Outlet 

Settings 

Flow Regime 

Subsonic 

Mass And Momentum 

Mass Flow Rate 

Mass Flow Rate 

MassFlow 

Mass Flow Rate Area 

As Specified 

Boundary - Bottom 

Type 

WALL 

Location 

Bottom 

Settings 

Heat Transfer 

Fixed Temperature 

Fixed Temperature 

9.3300e+02 [K] 

Mass And Momentum 

No Slip Wall 

Boundary - Sidel 

Type 

WALL 

Location 

Sidel 

Settings 

Heat Transfer 

Fixed Temperature 

Fixed Temperature 

9.3300e+02 [K] 

Mass And Momentum 

No Slip Wall 

Boundary - Side2 

Type 

WALL 

Location 

Side 2 

Settings 

Heat Transfer 

Fixed Temperature 

Fixed Temperature 

9.3300e+02 [K] 

Mass And Momentum 

No Slip Wall 
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Boundary - Top 

Type 

WALL 

Location 

Top 

Settings 

Heat Transfer 

Fixed Temperature 

Fixed Temperature 

9.3300e+02 [K] 

Mass And Momentum 

No Slip Wall 
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4. Solution Report 

Table 6. Boundary Flows for CFX 1 


Location 

Type 

Mass Flow 

Momentum 

X 

Y 

Z 

Bottom 

Boundary 

O.OOOOe+00 

4.3131e-10 

-2.7280e-01 

-6.8225e-04 

Inlet 

Boundary 

2.1000e-01 

-2.4106e-14 

3.7373e-14 

1.2185e-03 

Outlet 

Boundary 

-2.1000e-01 

1.3553e-07 

1.6927e-07 

1.5106e-03 

Sidel 

Boundary 

O.OOOOe+00 

-2.7280e-01 

3.0748e-10 

-6.8225e-04 

Side2 

Boundary 

O.OOOOe+00 

2.7280e-01 

3.3885e-10 

-6.8225e-04 

Top 

Boundary 

O.OOOOe+00 

4.1045e-10 

2.7280e-01 

-6.8225e-04 


Table 7. Forces and Torques for CFX 1 


Location 

Type 

X 

Y 

Z 

Bottom 

Pressure Force 

O.OOOOe+00 

2.7280e-01 

O.OOOOe+00 

Viscous Force 

-4.3132e-10 

5.663 le-07 

6.8225e-04 

Total Force 

-4.3132e-10 

2.7280e-01 

6.8225e-04 

Pressure Torque 

-7.2782e+01 

O.OOOOe+00 

-3.8641e-08 

Viscous Torque 

-1.0248e-03 

-1.6459e-07 

-6.4290e-10 

Total Torque 

-7.2783e+01 

-1.6459e-07 

-3.9284e-08 

Sidel 

Pressure Force 

2.7280e-01 

O.OOOOe+00 

O.OOOOe+00 

Viscous Force 

5.6630e-07 

-3.0748e-10 

6.8225e-04 

Total Force 

2.7280e-01 

-3.0748e-10 

6.8225e-04 

Pressure Torque 

O.OOOOe+00 

7.2782e+01 

5.6216e-08 

Viscous Torque 

1.1879e-07 

1.0248e-03 

3.9457e-10 

Total Torque 

1.1879e-07 

7.2783e+01 

5.661 le-08 

Side2 

Pressure Force 

-2.7280e-01 

O.OOOOe+00 

O.OOOOe+00 

Viscous Force 

-5.6761e-07 

-3.3885e-10 

6.8225e-04 

Total Force 

-2.7280e-01 

-3.3885e-10 

6.8225e-04 

Pressure Torque 

O.OOOOe+00 

-7.2782e+01 

-4.8739e-08 

Viscous Torque 

1.2992e-07 

-1.0253e-03 

-5.9973e-10 

Total Torque 

1.2992e-07 

-7.2783e+01 

-4.9338e-08 

Top 

Pressure Force 

O.OOOOe+00 

-2.7280e-01 

O.OOOOe+00 

Viscous Force 

-4.1046e-10 

-5.6760e-07 

6.8225e-04 

Total Force 

-4.1046e-10 

-2.7280e-01 

6.8225e-04 
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Pressure Torque 

7.2782e+01 

O.OOOOe+00 

3.9689e-08 

Viscous Torque 

1.0253e-03 

-1.6161e-07 

6.2895e-10 

Total Torque 

7.2783e+01 

-1.6161e-07 

4.0318e-08 
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5. Pictures 


Figure 1. Pressure Profile 
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Figure 2. Temperature Profile 
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Figure 3. Velocity Profile 
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APPENDIX D. ANSYS-CFX OPTIMIZATION REPORT FOR 50KW 

PLANT 


ANSYS 


Summary 


Project: 

SmgleDuct(Scaled) 

Date: 

5/23/2017 

Time: 

8:19:37 PM 

Product Version: 

Release 17.2 

Last Saved Version: 

Release 17.2 


Project Schematic View 


▼ A 


▼ B 


n ^ FUd Flow (CFX) 1 

2 ^ Geometry 

^ a 

2 Geometry 

^ a 

3 ^ Mesh 

^ a 

3 ^ Mesh 

^ a 

4 ^ Setup 


4 il Setup 

v' ^ 

5 Solution 

■y 

5 Solution 

v' ^ 

6 1^ Results 

y A 

6 ^ Results 


Quasi-20 Duct 


• 7 Cp5 Parameters 


3-D Single Duct (2 Duct System) 


(fa Parameter Set 


Response Surface OptimKation 


1 Iva^l _ 

2 13 Design of Experiments V 

3 ■ Response Surface V 

Optimization -y 

Response Surface Optimization 
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Files 


Name 

Size 

Type 

Date Modified 

Location 

i^SingleDuct(Scaled).x_t 

5 KB 

Geometry File 

12/2/2016 10:22:51 

AM 

F[:\MyDocs\Thesis\ 

Model 

_bpXl.cmdb 

6 MB 

CFX Mesh File 

2/20/2017 5:07:51 PM 

dp0\CFX\MECH 

(5)cFX.mshdb 

6 MB 

Mesh Database File 

2/20/2017 5:07:54 PM 

dp0\global\MECH 

A 

SingeDuct(Scaled).w 

bpj 

484 

KB 

Workbench Project 

File 

5/23/2017 12:24:36 

AM 

FI:\MyDocs\Thesis\ 

Model 


Design Points 


Name 

P5- 

MassFlow 

P2- 

DeltaT 

P3- 

Power 

P4- 

T3 

Retained 

Retained 

Data 

State 

Report 

Note 

Units 

kg s^-1 

K 

w 

K 





DPO 

(Curre 

nt) 

0.21 

282.91 

22653 

747.21 

yes 

✓ 



DP 1 

0.6184 

/ 

/ 





Created 

from 

Response 

Surface 


Outline of All Parameters 


ID 

Parameter Name 

Value 

Unit 

EIH Input Parameters 




E1H52)3-D Single Duct (2 Duct System) 




P5 

MassFlow 

0.21 

kg s''-l 

EIH Output Parameters 




EIBESIs-D Single Duct (2 Duct System) 




1^2 

DeltaT 

282.91 

K 
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P^P3 

Power 

22653 

W 

p3p4 

T3 

747.21 

K 


Response Surface Optimization 
Design of Experiments 

The Design of Experiments is the initial step building a Response Surface over the design 
space. This section describes the selected input parameters and their variation range, the 
chosen Design of Experiments type, and the generated Matrix of Experiments. 

Parameters 

The explored design space is defined by the range of variation of the following 1 input 
parameters. 


ID 

Name 

Classification 

Lower Bound 

Upper Bound 

P5 

MassFlow 

Continuous 

0.01 [kgs^-1] 

1 [kg s^-l] 


Design of Experiments Properties 


Property 

Value 

Design of Experiments Type 

Central Composite Design 

Design Type 

Auto Defined 


Matrix of Experiments 

All of the 5 points are up-to-date. 

The full Matrix of Experiments is provided in the Appendices. 
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Parameters Parallel Chart 


Parameters Parallel Chart 

407 33 2 9437 


933.01 



0 0099999 
P5 (kgs^-11 


117.26 
P2 [K1 


0 18507 
P3 (xlO*) (W) 


495 03 
P4 (K1 
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Design Points vs Parameter Chart 


Design Points vs Parameter mm 



Points Name 


Response Surface 

The Response Surface is a meta-model built from the Design of Experiments for an 
efficient exploration of the design space. This section describes the selected type of meta¬ 
model, including its properties, the obtained quality, and the generated Response Points 
and charts. 

Response Surface Properties 


Property 

Value 

Response Surface Type 

Genetic Aggregation 

Random Generator Seed 

0 

Maximum Number of Generations 

12 


89 
















Display Level 

Final 

Generate Verification Points 

No 


Goodness of Fit 


The Goodness of Fit report provides charts and metrics to understand how each output 
parameter is approximated by the response surface. 

Details of ‘Goodness Of Fit’ 

Predicted vs Observed Scatter Chart 


Predicted vs Observed - Normalized Values msm 




P2 - DeltaT 

P3 - Power 

P4-T3 

Coefficient of Determination (Best Value = 1) 





if 

it 

ir 

Learning Points 
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Cross-Validation on Learning Points 

'A’l 


irki 

Root Mean Square Error (Best Value = 0) 




Learning Points 

7.8215E-05 

4.4538E-05 

7.1423E-06 

Cross-Validation on Learning Points 

0.00037587 

0.00033711 

L3094E-05 

Relative Maximum Absolute Error (Best Value = 0%) 




Learning Points 

Ao 

Ao 

Ao 

Cross-Validation on Learning Points 

'*'*’0.00066882 

'*'*6.1964E-06 

'**’L0682E-05 

Relative Average Absolute Error (Best Value = 0%) 




Learning Points 

Ao 

Ao 

Ao 

Cross-Validation on Learning Points 

*"*0.00022649 

*"*L7665E-06 

*"*6.3594E-06 


Minimum and maximum values 

This section reports the minimum and maximum values for each output parameter. These 
values are approximations found by the Min-Max Search on the Response Surface. 


Name 

Minimum value 

Maximum value 

P2 - DeltaT (K) 

117.26 

407.33 

P3 - Power (W) 

1850.7 

29731 

P4 - T3 (K) 

495.04 

933 


The input parameter values corresponding to the minimum and maximum values are 

provided in the Appendices. 

Response Points 

The Response Points provide the output parameter values obtained by evaluating the 
Response Surface. This section lists all the Response Points and their associated charts. 


Name 

P5 - MassFlow (kg 

P2 - DeltaT 

(K) 

P3 - Power 
(W) 

P4-T3 

(K) 

Report 
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Response 

Point 

0.505 

187.07 

29424 

604.63 



Details of ‘Response Point’ 
Response Chart 


Response Chart for P3 - Power msfs 
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Local Sensitivity Chart 



P2 - DeltaT P3 - Power P4 - T3 

Output Parameters 
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Local Sensitivity Curves Chart 

Axes Range = Use Min Max of the Output Parameter 


Local Sensitivity Curves P5 - MassFlow vs P2 - Delta! /MSIS 



94 





Spider Chart 



Optimization 

The Optimization is based on Response Surface evaluations. This section describes the 
chosen Optimization type and the generated candidates and charts. 

Parameters 

The explored design space is defined by the range of variation of the following 1 input 
parameters. 


ID 

Name 

Classification 

Lower Bound 

Upper Bound 

P5 

MassFlow 

Continuous 

0.01 [kgs^-1] 

1 [kgs^-l] 
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Properties 


Property 

Value 

Method Name 

Screening 

Verify Candidate Points 

No 

Number of Samples 

1000 

Maximum Number of Candidates 

3 


Definition 

Objectives and constraints of the optimization study 


P5 - MassFlow (kg s'^-l) 

P3 - Power (W) 

Domain 

[0.01;1] 


Objective 


Maximize 


Description 


Optimization 

Method 

Screening 

The Screening optimization method uses a simple approach based on sampling and sorting. It 
supports multiple objectives and constraints as well as all types of input parameters. Usually it 
is used for preliminary design, which may lead you to apply other methods for more refined 
optimization results. 

Configuration 

Generate 1000 samples and find 3 candidates. 

Status 

1 _ 

Converged after 1000 evaluations. 


Status 


Property 

Value 

Converged 

Yes 

Number of Evaluations 

1000 

Number of Failures 

0 

Size of Generated Sample Set 

1003 
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Number of Candidates 3 


Candidates 

Candidates of the optimization study 



P5 - MassFIow (kg 

P2 - DeltaT 

(K) 

P3 - Power 
(W) 

P4-T3 

(K) 

Report 

Candidate Point 

1 

0.6176 

165.24 

'*’*29731 

570.75 


Candidate Point 

2 

0.67875 

155.31 

**29664 

555.37 


Candidate Point 

3 

0.74013 

146.41 

**29484 

541.51 



The table of Rating Values is provided in the Appendices. 
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History Chart 
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History Chart 
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Tradeoff Chart 


Number of Pareto Fronts to Show = 995/995 


TradeaffIcHart P? - Delta? vs P3 ■ Power 
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Samples Chart 


Mode = Pareto Fronts 

Coloring method = per Pareto Front 

Number of Pareto Fronts to Show = 995/995 
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Samples Chart 
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Sensitivities Chart 
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Candidates Chart 


Coloring method = per Candidate 


Candidate Points 
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Appendix 

Matrix of Experiments (Response Surface Optimization system) 


Name 

P5 - MassFlow (kg s'^-l) 

P2 - DeltaT (K) 

P3 - Power (W) 

P4 - T3 (K) 

Report 

1 

0.505 

187.07 

29424 

604.63 


2 

0.01 

407.33 

1850.7 

933 


3 

1 

117.26 

27959 

495.04 


4 

0.2575 

260.01 

24557 

713.76 


5 

0.7525 

144.73 

29437 

538.88 
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Output Parameter Minimums (Response Surface Optimization system) 


Name 

P5 - MassFlow (kg s'^-l) 

P2 - DeltaT (K) 

P3 - Power (W) 

P4 - T3 (K) 

P2 - DeltaT 

1 

117.26 

27959 

495.04 

P3 - Power 

0.01 

407.33 

1850.7 

933 

P4-T3 

1 

117.26 

27959 

495.04 


Output Parameter Maximums (Response Surface Optimization system) 


Name 

P5 - MassFlow (kg s'^-l) 

P2 - DeltaT (K) 

P3 - Power (W) 

P4 - T3 (K) 

P2 - DeltaT 

0.01 

407.33 

1850.7 

933 

P3 - Power 

0.6176 

165.24 

29731 

570.75 

P4-T3 

0.01 

407.33 

1850.7 

933 


Table of Rating Values (Response Surface Optimization system) 


P3 - Power (W) 

Objective 

Maximize 

A 

[27407; 29731 ] 

★★ 

[ 22761 ; 27407 ] 

★ 

[ 18114 ; 22761 ] 

— 

[ 13467 ; 18114] 

X 

[ 8820.7 ; 13467 ] 

XX 

[4174 ; 8820.7 ] 

X 

XX 

[ 1850.7 ; 4174] 
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APPENDIX E. ANSYS-CEX OPTIMIZATION RAW DATA EOR 

50KW PLANT 


Mass Flow [kg/ 
s] 

0.02 

0.1025 

0.185 

0.2675 

0.35 

0.4325 

0.515 

0.5975 

0.68 

0.7625 

0.845 

0.9275 

1.01 

1.0925 

1.175 

1.2575 

1.34 

1.4225 

1.505 

1.5875 

1.67 

1.7525 

1.835 

1.9175 

2 


T3 [K] 

933 

898.6047457 

862.0005805 

824.1123246 

785.8689574 

748.3235337 

713.7579984 

691.6929407 

671.6875067 

652.98166 

635.536148 

619.3964533 

604.6250075 

591.2289457 

579.0219231 

567.8161804 

557.4652728 

547.8526508 

538.8839875 

530.4818445 

522.5818596 

515.1299579 

508.0802721 

501.3935634 

495.0360063 


dT[K] 

407.3290384 

377.5018905 

349.7800669 

324.113039 

300.4774403 

278.9230139 

260.0088857 

244.6514758 

230.9024898 

218.3910802 

206.9727605 

196.5559324 

187.0721207 

178.4503371 

170.574234 

163.3392515 

156.6625267 

150.4770975 

144.7279626 

139.3693151 

134.3625463 

129.6747742 

125.2777347 

121.1469344 

117.2609916 


Power [kW] 

3.701339955 

17.23358692 

27.90234773 

35.49245444 

41.08261961 

45.53495101 

49.11380011 

51.97397773 

54.22868854 

55.96955214 

57.27484372 

58.21355049 

58.84759985 

59.23398199 

59.42401775 

59.45684005 

59.3618072 

59.1614127 

58.87320006 

58.51106224 

58.08614696 

57.60750197 

57.08254418 

56.51740676 

55.91720001 
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